A semi- schematic model for the center of mass dynamics in 
supercooled molecular liquids 



Linda Fabbian, Francesco Sciortino, Filippo Thiery and Piero Tartaglia 
Dipartimento di Fisica and Istituto Nazionale per la Fisica della Materia, Universitd di Roma La 

Sapienza, P.le Aldo Moro 2, 1-00185, Roma, Italy 
(February 7, 2008) 



Abstract 

We introduce a semi-sciiematic mode-coupling model to describe the slow 
dynamics in molecular liquids, retaining explicitly only the description of the 
center of mass degrees of freedom. Angular degrees of freedom are condensed 
in a g-vector independent coupling parameter. We compare the time and q- 
dependence of the density fluctuation correlators with numerical data from a 
250 ns long molecular dynamics simulation. Notwithstanding the choice of a 
network-forming liquid as a model for comparing theory and simulation, the 
model describes the main static and dynamic features of the relaxation in a 
broad g-vector range. 
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The slow dynamics in supercooled liquids has recently become a field of increasing in- 
terest [Q and a large amount of knowledge has been accumulated from theoretical work , 
experiments |^,^ and simulations It has now become apparent that the ideal mode- 

coupling theory (MCT) for structural relaxation contains the essential ingredients for 
the description of the molecular dynamics in supercooled states, down to a cross-over tem- 
perature Tc, where "hopping" processes become dominant [pjp^. 

Although the MCT was originally developed for simple (atomic) liquids, it appears to be 
able to rationalize in a coherent way experimental 0,|ll],|12[ and numerical simulation results 
,T3,n] also for molecular liquids, i.e. for molecules with anisotropic interaction potentials. 



The lack of theoretical predictions for the behavior of the center of mass (COM) and angu- 
lar degrees of freedom in molecular liquids has been compensated by the use of schematic 
models, where the spectrum of different g-vectors and orientations is condensed into one or 
two representative correlators. The exact solution of these schematic models and the com- 
parison with experiments [0,11^ and simulations [T^ provide compelling evidence that even 
molecular liquids can be described in a MCT framework. A full and detailed study assessing 
the ability of MCT to predict not only the general features of supercooled molecular liquids, 
but also the g-dependence of the COM correlators is highly desirable. 
The possibility of describing theoretically the slow dynamics in supercooled molecular liq- 
uids is of fundamental importance for assessing the limitations of the ideal MCT, in making 
contact with experiments, and in tackling the problem of the relation between translational 
and rotational motions in supercooled states. Theoretical work has recently been under- 
taken |lT7|-pO[] in this direction and appears to offer in future the possibility of a detailed 
comparison between experiments and theory. 

In this Letter we propose a semi- schematic mode-coupling model to describe the COM 
dynamics in molecular liquids and we compare the predictions of the model with the COM 
collective correlation functions, calculated from a long molecular dynamics (MD) simulation 
for a model of a network forming molecular liquid |]T1|. For the first time, the comparison is 



performed at the highest possible level of detail, i.e. comparing not only static quantities, 
like the g-dependence of the non-ergodicity factors, but also the time dependence of the 
correlators in the entire q range. To our knowledge, such level of detailed comparison has 
never been achieved, not even for the simple liquid case. The model we propose is a modifi- 
cation of the MCT equations for simple hquids. It takes into account the coupling between 
rotational and translational degrees of freedom in a phenomenological way, by introducing 
a g-independent coupling parameter xr which measures the increase in the strength of the 
memory function m,q{t) induced by the roto-translational coupling. For = 1 we recover 
the usual MCT equations for simple liquids, while Xi? > 1 models a system in which the ro- 
tational degrees of freedom slow down the COM dynamics. In contrast to schematic models 



Tql we retain all the g-dependence information for the COM motion In our MCT 



model the COM dynamics in the a-relaxation region is entirely controlled by the value of the 
COM structure factor Sg, by the number density n and by the value of xr- The latter can 
be considered as a fitting parameter which controls the ideal glass transition temperature. 
The only T or P dependence in the dynamics arises from the T and P dependence of Sq 
and n. 

With the choice for m,q{t) discussed above, the relaxation of the density- density correlation 
functions 4>g{t) in the a-region is described, in terms of a rescaled time t = t/T [^, by 
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a system of coupled integro-differential equations: 
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Eq.|l| and Eq.^ are the usual MCT expressions while Eq.|^ differs from the one for simple 
liquids only by the presence of the g- independent multiplicative factor xr- the early a- 
relaxation region the asymptotic analytic solution of Eq.|l] follows, to leading order, a power 
law (the von Schweidler law) 
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where fg, the so-called non-ergodicity factor, can be calculated solving the coupled integral 
equations 



mq[fk]. 



The exponent b in Eq. ^ is related to the exponent parameter A by 

r(i + 26) 

where F is the Euler gamma function. A is defined by 
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where e'^ and e'^ are the right and left eigenvectors of the maximum eigenvalue of the stability 
matrix 
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evaluated at the ideal transition temperature Tc. Details on the calculation of the critical 
amplitudes h^^ and h^^^ can be found in Ref. |^ . 

In the late a-region it is not possible to solve analytically Eq.|l|, but the numerical solution 
is well fitted by a stretched exponential law (Kohlrausch- William- Watts law) 
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where < 1. 

We compare the predictions of the semi-schematic model (Eqs.|I|-0) with the corresponding 
correlators evaluated in a long MD simulation [Q. The simulation was performed us- 
ing the SPC/E potential [^], a model developed to simulate the behavior of liquid water. 
It treats the molecule as a rigid unit and models the pair interactions as a sum of elec- 
trostatic and Lennard Jones terms. Previous studies have shown that SPC/E reproduces 
semi-quantitatively the static and dynamic properties of liquid water p5[. In water, an 
open-structured liquid with an average coordination number around 4, the slowing down of 
the dynamics on supercooling is not related to packing constraints but it is driven by the 
formation of a tetrahedral network of strong and highly directional bonds. The role of the 
cages in reducing the mobility of the molecules is played by the hydrogen bonds which freeze 
the rotational motion of the molecules and therefore hinder their translational motion. The 
steric cages which slow down the dynamics in simple liquids are replaced by "energetic" 
cages. 

It has been shown numerically that the SPC/E dynamics in supercooled states is character- 
ized by a two step relaxation process, by a self-similar region in the early a-relaxation region 
(with a von-Schweidler exponent b = 0.50 ± 0.05) and by a stretched exponential decay in 
the long time a regime. The relaxation times display a power law divergence at the critical 
temperature = (200 ± 2)K, which has been identified with the ideal MCT transition 



temperature |T^. Moreover, the rotational correlators follow closely the COM correlators, 
supporting the view that hydrogen bonds form strong energetic cages and that molecules 
can not rotate within a cage. The process of breaking and reforming cages, characteristic of 
the a-relaxation region is the bottleneck for both COM and rotational diffusion. 
We solve numerically Eqs.|I|-^ on a grid of 300 equispaced q vectors, using as input the COM 
number density {n = 32.3 nm~^) and the COM Sg evaluated from the MD data. We study 
the model for different values of the coupling xr- When xr = 1-93, the ideal glass transition 
temperature in the model coincides with the Tc detected in SPC/E water [|I4|. Since with 
this choice of xr the ideal glass transition temperatures in the model and in SPC/E are 
the same, we can compare, without using further fitting parameters, results from the MD 
simulations and theoretical predictions for all g-vectors. 

The comparison for the g-dependence of fq, h^^^ and h^^^ in Eq.^ is shown in Fig. [I| 



p7| . For reference, we also show 5*^, the only g-dependent input in the theory. As in the 
simple liquid case, the static quantities show oscillations in phase with Sq. We stress, in 
passing, the peculiar shape of the COM structure factor Sq for SPC/E water when compared 
to the Sq of simple liquids. The presence of the pre-peak around q = 18nm~^ (the analog 
of the so-called first sharp diffraction peak in silica) is the hallmark of medium range order 
characteristic of network- forming liquids [p^ , ^] . 

The agreement between numerical and theoretical values for fq and h^q^ is striking although 
deviations at large q are clearly seen; the difference between theory and simulations for large 
values of q {i.e. small distances) is not surprising since this is the region where the system 
is more sensitive to the detailed geometry of the molecules, completely neglected by the 
model. The agreement for h^^^ is not as good. We recall that the numerical values for h^^^ 
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were obtained by a fitting of (pq{t) to a second order polynomial in (t/r)^ (Eq.§) and thus 
the fitted values are affected by higher order corrections; on the contrary the theoretical h^^^ 
is the exact coefficient of the second order term (t/r)^''. 

From the calculated fg values, we evaluate the exponent parameter A by numerical integra- 
tion of Eq. 0. We obtain A = 0.79 which yields, through Eq. |^, 6 = 0.49. The value for 
the exponent b, which controls the self-similar dynamics in the early a-region, is in excellent 
agreement with the one calculated by MD, i.e. b = 0.50 it 0.05. Thus, the discrepancies 
between theory and numerical data at large g-values do not affect the value of the von 
Schweidler exponent b. Indeed, the leading contributions to the integral in Eq.[7| come from 
the maxima of 5*^, i.e. from the values at the pre-peak and at the main peak, where our 
model provides the correct results. 

To extend the comparison to the late a-region, we solve Eq|^ for all 300 g-vectors and 
we fit, for each q, the time evolution of the MD and theoretical correlators to a Kohlrausch 
stretched exponential according to Eq.^ The g-dependence of the fit parameters (amplitude 
, relaxation time and stretching exponent f3^) is compared with the analogous quan- 
tities obtained from the MD data in Fig.^ Again, we find a remarkable agreement between 
the theoretical predictions and the MD data in the low q region, both in the amplitudes and 
in the relaxation times. Instead, the parameters f3^ show oscillations in phase with 5*^ but 
the theoretical and MD values differ up to 25%. This indicates that in the MD data, the 
a-relaxation dynamics covers a larger time range. The differences in (3^ can be tentatively 
associated with the drastic reduction in the number of coupled correlators retained in the 
present model, which neglects the angular correlators. 

In summary, a simple modification of the MCT equations for simple liquids is able to describe 
the a-relaxation process in SPC/E water, a model for a molecular network- forming liquid 
with strong directional bonds. By using as input the fixed value of n and the COM structure 
factor 5*^, the theory is able to describe both the q- and t-dependence of the density-density 
correlations in a satisfactory way. We stress again that except for the value of xr and for 
a g-independent scaling of time relating t with real time |2^, no fitting parameters are in- 



volved in the comparison. The model can be extended to describe the self-particle dynamics 



and the collective transverse currents \31 



The approach presented here is a first step towards a full MCT description, which will nec- 
essarily include orientational correlations, i.e. the geometry of the molecule [p!7|-[20|]. Our 



semi-schematic approach provides evidence that the slow dynamics in complex molecular 
liquids can be quantitatively described by a very simple modification of the MCT equations 
for simple liquids, thus retaining the simplicity of the COM description. 
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FIGURES 



FIG. 1. Comparison between theory (solid lines) and MD simulations (symbols) for the 



non-ergodicity parameter fq and the critical amplitudes h 



(1) 



and h'^q^ 



The center of mass structure factor Sq is shown as a reference. MD data are from Ref. 



in Eq.|. (See Ref. 

0- 



FIG. 2. Fitting parameters to the Kohlrausch- William- Watts law (Eqj9|) as obtained by the 
solution of Eq.|^ (solid lines) and by the MD data (symbols). The time scale for the MD values of 

is 1 ps. In the theoretical values of the time scale has been chosen in an appropriate way 
as explained in Ref. [27|. MD data are from Ref. |14|. 
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FIG. 1. L. Fabbian et al. 
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FIG. 2. L. Fabbian et al. 



